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_ ■ ABSTRACT 

. High-energy collisions of cosmic-ray nuclei with interstellar gas are believed to be the mechanism 

' producing the majority of cosmic ray antiprotons. Due to the kinematics of the process they are created 

with a nonzero momentum; the characteristic spectral shape with a maximum at ~ 2 GeV and a sharp 
decrease towards lower energies makes antiprotons a unique probe of models for particle propagation 
in the Galaxy and modulation in the heliosphere. On the other hand, accurate calculation of the 
£N1 ■ secondary antiproton flux provides a "background" for searches for exotic signals from the annihilation of 

supersymmetric particles and primordial black hole evaporation. Recently new data with large statistics 
on both low and high energy antiproton fluxes have become available which allow such tests to be 
performed. We use our propagation code GALPROP to calculate interstellar cosmic-ray propagation 
for a variety of models. We show that there is no simple model capable of accurately describing the 
■ whole variety of data: boron/carbon and sub-iron/iron ratios, spectra of protons, helium, antiprotons, 

^—1 | positrons, electrons, and diffuse 7-rays. We find that only a model with a break in the diffusion coefficient 

plus convection can reproduce measurements of cosmic-ray species, and the reproduction of primaries 
(p, He) can be further improved by introducing a break in the primary injection spectra. For our best-fit 
model we make predictions of proton and antiproton fluxes near the Earth for different modulation levels 
1 ' and magnetic polarity using a steady-state drift model of propagation in the heliosphere. 

2 ' Subject headings: diffusion — convection — elementary particles — nuclear reactions, nucleosynthesis, 

I abundances — cosmic rays — ISM: general — Galaxy: general — cosmology: theory 

C$ • — dark matter 

• ^h. 1. INTRODUCTION cosmological abundance and be present in our own Galaxy. 

' T. x , r , , *~,-pj , . , , , , -, T-, , , A pair of stable WIMPs can annihilate into known parti- 

Most of the CR antiprotons observed near the Earth are , . , . , . . .. ., ... . - tti — TTV /r~i — 1 • 

53 , , . j j • ,,. . c ,. ^r, cles and antiparticles making it possible to infer WIMPs m 

secondaries produced m collisions of energetic CR parti- , • , , , , , r , • •, ■ 

cles with interstellar gas (e.g., Mitchell et al. 1996). Due to the Galactic halo by the products of their annihilations 

, . , , c '. , PBHs may have formed m the early Universe via initial 

the kinematics ot this process, the spectrum ot antiprotons J . . . ■' 

, , . . • , • .. r ,1 • density fluctuations, phase transitions, or the collapse of 

has a unique shape distinguishing it from other cosmic-ray . . , TT . " , _, , , , . , ,. . „ 

species. It peaks at about 2 GeV decreasing sharply to- T^V^^T^ ' 1 J ' 

, , . T ,,,,, , ,. , Onto 1996 . Black holes can emit particles and evaporate 

wards lower energies. In addition to secondary antiprotons rr mi • • • ^ 

c i • i n due to quantum effects. I he emission rate is generally too 

there are possible sources of primary antiprotons; those . , . . ° J , , 

i pi j ■ j ,i j 1 ,-i -i • low to be observable, but it increases as the black hole 

most often discussed arc the dark matter particle annini- , ' . , . „„ TT . . 

, , , , . r j- i ii i i i ^nnTT \ mass decreases. Ihe only observable PBHs are those that 

lation and evaporation of primordial black holes (PBHs). . „ , . . . . . 

m, , , . r ,i i i , , have a mass small enough to produce a burst of particles 

1 he nature and properties ot the dark matter that con- or- r 

, ., , . .„ c . c ,, „ , j j . as they evaporate. 

stitute a sigmhcant fraction of the mass of the Universe T J ^ , J 1 , 1 , . t . 

i ii r ,i i j / m • i i In recent years, new data with large statistics on both 

have puzzled scientists for more than a decade ( Irimble . . / ' n ° ., 

i non a t, i nnn\ a j-t, c a a i l° w a n0 - high energy antiproton huxes have become avail- 

1989; Ashman 1992). Among the favored dark matter can- able fHof ct al 1996 B s et 1 1999 O 't t 1 2000 

didates are so-called weakly-interacting massive particles _ \ .. '. „„'„ „ } m ' . „'„ l rl r , , .' ,' 

/™mv,m \ v. ■ i c ii c a • Bergstrom et al. 2000; Maeno et al. 2001; Stochai et al. 

(WIMPs), whose existence follows from supersymmetric » , , „™ n , „ i , 

j t ic ■ t v ■ i i • e ■ 2001; Asaoka et al. 2001a) that allow us to test models 

models (lor a review see Jungman, Kamionkowski, & Gn- ' . , , • A i 

est 1996). Such particles, if stable, could have a significant of CR P ro P a S atl0n and hehosphcric modulation. A probe 
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Table 1 

Summary of recent CR antiproton calculations. 



Model 


Antiproton 
production 


Ambient CR 
spectrum 


Propagation 
model 


Reacceleration / 
Convection 


Modulation 
in heliosphere 


Tertiary 
antiprotons 


Bergstrom, Edsjo, & Ullio (1999) 


Tan & Ng 


scaled LIS 


simplified 
diffusion 


n/y 


force-field 


y 


Bieber et al. (1999) 


DTUNUC 
event generator 


LIS, in the 
whole Galaxy 


leaky-box 


n/n 


steady-state 
drift modcl a 


y 


oottmo et al. (1990) 


lan 6i JNg 


Llo, in the 
whole Galaxy 


two zone 
diffusion 


n/n 


lorcc- field 


n 


Donate et al. (2001) 


Tan & Ng, 
DTUNUC 


LIS, in the 
whole Galaxy 


two zone 
diffusion 


y/y 


force-field 


y 


Simon, Molnar, & Rocslcr (1998) 


DTUNUC 
event generator 


LIS, in the 
whole Galaxy 


leaky-box 


n/n 


force-field 


ii 


Present work 


pp: Tan & Ng 
pA: DTUNUC 


propagated 


3D spatial 
grid, realistic 


y/y 


steady-state 
drift model 


y 



diffusion 

a Parameters are fixed from the assumed LIS proton spectrum. 



to measure low energy particles in interstellar space may 
also become reality in the near future (Mewaldt & Liewer 
2001). Additionally, accurate calculation of the secondary 
antiproton flux provides a "background" for searches for 
exotic signals such as WIMP annihilation or PBH evapo- 
ration. 

Despite numerous efforts and overall agreement on the 
secondary nature of the majority of CR antiprotons, pub- 
lished estimates of the expected flux significantly differ 
(see, e.g., Fig. 3 in Orito et al. 2000). Calculation of the 
secondary antiproton flux is a complicated task. The ma- 
jor sources of uncertainties are three fold: (i) incomplete 
knowledge of cross sections for antiproton production, an- 
nihilation, and scattering, (ii) parameters and models of 
particle propagation in the Galaxy, and (iii) modulation 
in the heliosphere. While the interstellar antiproton flux 
is affected only by uncertainties in the cross sections and 
propagation models, the final comparison with experiment 
can only be made after correcting for the solar modulation. 
Besides, the spectra of CR nucleons have been directly 
measured only inside the heliosphere while we need to 
know the spectrum outside, in interstellar space, to com- 
pute the antiproton production rate correctly. The basic 
features of the recent models are summarized in Tabic 1. 

For the antiproton production cross section two options 
exist: a semiphenomenological fit to the data by Tan & 
Ng (1983a, b) which has been used for almost two decades 
because of lack of new data and the Monte Carlo event 
generator DTUNUC (Ferrari et al. 1996; Roesler, Engel, & 
Ranft 1998). Both the parametrization and the DTUNUC 
code describe well the available data on antiproton pro- 
duction in pp-collisions. For interactions involving heav- 
ier nuclei there arc no measurements so far. However, 
the DTUNUC model provides a reasonable description of 
hadron-nucleus and nucleus-nucleus interactions and thus 
it can be used for estimates of antiproton production in 
proton- nucleus collisions. Additionally, "tertiary" antipro- 
tons, inelastically scattered secondaries, are the main com- 



ponent at low energies. This component has only been 
included in the most recent models. 

To calculate particle propagation one must choose be- 
tween the easy-to-apply but non-physical "leaky-box" 
model, and a variety of diffusion models. The leaky-box 
model reduces the problem to derivation of the path length 
distribution (which in turn remains the major source of un- 
certainty as discussed, e.g., by Simon, Molnar, & Roesler 
1998). Diffusion models are theoretically more physical, 
but also differ by the degree to which they reflect reality. 
The ultimate goal is to develop a model which is consistent 
with many different kinds of data available: measurements 
of CR species, 7-rays, and synchrotron emission. 

Heliosphcric modulation is not yet understood in detail. 
The first and yet popular force-field approximation (Glee- 
son & Axford 1968) does not work at low energies. More 
sophisticated drift models reflect our still incomplete cur- 
rent knowledge. The main problem here is that the modu- 
lation parameters are determined based on the assumed ad 
hoc interstellar nucleon spectrum. The antiproton/proton 
ratio, which is often calculated, is even more uncertain at 
low energies and may vary by an order of magnitude over 
the solar cycle (Labrador & Mewaldt 1997). It is therefore 
the antiproton spectrum itself which should be compared 
with calculations. On the other hand, a reliable calcu- 
lation of the antiproton spectrum would allow the study 
of charge sign dependent effects in the heliosphere. That 
could, in turn, help to derive the low-energy part of the 
local interstellar (LIS) spectrum of nucleons and help to 
establish the heliosphcric diffusion coefficients more accu- 
rately. 

We have developed a numerical method and correspond- 
ing computer code GALPROP 3 for the calculation of 
Galactic CR propagation in 3D (Strong & Moskalenko 
1998). The code has been shown to reproduce simulta- 
neously observational data of many kinds related to CR 
origin and propagation (for a review see Strong & Mos- 
kalenko 1999; Moskalenko & Strong 2000). The code has 



3 Our model including software and datasets is available at http://www.gamma.mpe-garching.mpg.de/~aws/aws.html 
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been validated on direct measurements of nuclei, antipro- 
tons, electrons, and positrons, and astronomical measure- 
ments of 7-rays and synchrotron radiation. These data 
provide many independent constraints on model parame- 
ters. 

The code is sufficiently general that new physical effects 
can be introduced as required. The basic spatial propa- 
gation mechanisms arc diffusion and convection, while in 
momentum space energy loss and diffusive reacceleration 
are treated. Fragmentation, secondary particle produc- 
tion, and energy losses are computed using realistic distri- 
butions for the interstellar gas and radiation fields. Re- 
cently the code has been entirely rewritten in CH — h with 
many new features added (Strong & Moskalenko 2001). 

In this paper we use the GALPROP code for accurate 
calculation of production and propagation of secondary 
antiprotons and positrons. We explore the dependence of 
the antiproton and positron fluxes on the nucleon injection 
spectrum and propagation parameters. The antiproton 
production cross section is calculated using the pp pro- 
duction cross section (Tan & Ng 1983b) and DTUNUC 
nuclear factors from Simon ct al. (1998) to calculate p pro- 
duction in pA collisions or the pp production cross section 
scaled appropriately with beam/target atomic numbers. 
Inelastic scattering to produce "tertiary" antiprotons and 
"secondary" protons is taken into account. The calculated 
proton and antiproton LIS and p/p ratio are modulated us- 
ing the steady-state drift model. Our preliminary results 
have been reported in Moskalenko et al. (2001b). 

2. BASIC FEATURES OF THE GALPROP MODELS 

The cylindrically symmetric GALPROP models have 
been described in detail elsewhere (Strong & Moskalenko 
1998); here we summarize their basic features. 

The models are three dimensional with cylindrical sym- 
metry in the Galaxy, and the basic coordinates are (R, z,p) 
where R is Galactocentric radius, z is the distance from the 
Galactic plane and p is the total particle momentum. In 
the models the propagation region is bounded by R — Rh 1 
z = ±z h beyond which free escape is assumed. 

The propagation equation we use for all CR species is 
written in the form: 

dih c) Q ' 

-f t = q(r,p) + V • (D xx Vi/> W) + Q- p P 2 D PP Q- p j 2 V 



d_ r 

dp 



pip 



f(V-V)V 



—if) - —if) , 

Tf T r 



(1) 



where ip = ip(r,p,t) is the density per unit of total par- 
ticle momentum, ip(p)dp = 47rp 2 /(p) in terms of phase- 
space density /(p), q(r,p) is the source term, D xx is the 
spatial diffusion coefficient, V is the convection velocity, 
reacceleration is described as diffusion in momentum space 
and is determined by the coefficient D pp , p = dp/dt is 
the momentum loss rate, 77 is the time scale for frag- 
mentation, and T r is the time scale for radioactive de- 
cay. The numerical solution of the transport equation is 
based on a Crank-Nicholson (Press et al. 1992) implicit 
second-order scheme. The three spatial boundary condi- 
tions ip(Rh, z,p) — if)(R, ±Zh,p) — are imposed on each 
iteration, where we take Rh — 30 kpc. 



For a given Zh the diffusion coefficient function 
of momentum and the reacceleration or convection pa- 
rameters is determined by boron-to-carbon (B/C) ratio 
data. The spatial diffusion coefficient is taken as D xx — 
f3Do(p/ po) S if necessary with a break (5 — Si below rigid- 
ity po, S = 62 above rigidity po). The injection spectrum 
of nucleons is assumed to be a power law in momentum, 
dq(p)/dp oc for the injected particle density. 

Reacceleration provides a mechanism to reproduce the 
B/C ratio without an ad- hoc form for the diffusion co- 
efficient. Our reacceleration treatment assumes a Kol- 
mogorov spectrum with 5 = 1/3 or value close to this. 
For the case of reacceleration the momentum-space dif- 
fusion coefficient D pp is related to the spatial coefficient 
D xx (Berezinskii et al. 1990; Seo & Ptuskin 1994) via the 
Alfven speed va- 

The convection velocity (in z-direction only) V(z) is as- 
sumed to increase linearly with distance from the plane 
(dV/dz > for all z). This implies a constant adia- 
batic energy loss; the possibility of adiabatic energy gain 
(dV/dz < 0) is not considered. The linear form for V(z) is 
consistent with cosmic-ray driven MHD wind models (e.g., 
Zirakashvili ct al. 1996). The velocity at z = is a model 
parameter, but we consider here only V(0) = 0. (V(0) = 
is physically plausible: it is required to avoid a disconti- 
nuity, since by symmetry a non-zero velocity would have 
to change sign at z = 0.) 
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1. — Number density distributions of 2xH2 (solid), H I 



(dashes), and H II (dots) in the Galaxy. Shown are the plots for 
z = 0,0.1,0.2 kpc (decreasing density). Number density of H2 at 
2 = 0.2 kpc from the plane is very low and is not shown in the plot. 



The interstellar hydrogen distribution uses H 1 and CO 
surveys and information on the ionized component; the he- 
lium fraction of the gas is taken as 0.11 by number. The 
H2 and H 1 gas number densities in the Galactic plane 
are defined in the form of tables, which are interpolated 
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linearly. The extension of the gas distribution to an ar- 
bitrary height above the plane is made using analytical 
approximations. The distributions of 2xH 2 , H I, and H II 
are plotted on Figure 1 for z = 0,0.1,0.2 kpc. The code 
uses the densities averaged over the z-grid using a 0.01 
kpc step. More details about the gas distribution model 
are given in Appendix A. 

The distribution of cosmic-ray sources is chosen to re- 
produce the cosmic-ray distribution determined by anal- 
ysis of EGRET 7-ray data (Strong & Mattox 1996) and 
was described in Strong & Moskalenko (1998). 

Energy losses for nucleons by ionization and Coulomb 
interactions are included, and for electrons by ionization, 
Coulomb interactions, bremsstrahlung, inverse Compton, 
and synchrotron. 

Positrons and electrons (including secondary electrons) 
are propagated in the same model. Positron production 
is computed as described in Moskalenko & Strong (1998), 
that paper includes a critical reevaluation of the secondary 
7r ± - and J^ ± -meson decay calculations. 

Gas-related 7-ray intensities are computed from the 
emissivities as a function of (R, z,E 7 ) using the column 
densities of H 1 and H 2 . The interstellar radiation field 
(ISRF), used for calculation of the inverse Compton (IC) 
emission and electron energy losses, is calculated based 
on stellar population models and COBE results, plus the 
cosmic microwave background. 

An overview of our previous analyses is presented in 
Strong & Moskalenko (1999) and Moskalenko & Strong 
(2000) and full results for protons, helium, positrons, and 
electrons in Moskalenko & Strong (1998) and Strong, Mos- 
kalenko, & Reimcr (2000). The evaluations of the B/C and 
10 Be/ 9 Be ratios, diffusion/convection and reacceleration 
models and full details of the numerical method are given 
in Strong & Moskalenko (1998). Antiprotons have been 
previously discussed in the context of the "hard interstel- 
lar nucleon spectrum" hypothesis in Moskalenko, Strong, 
& Reimcr (1998) and Strong et al. (2000). Our results 
for diffuse continuum 7-rays, synchrotron radiation, and a 
new evaluation of the ISRF are described in Strong et al. 
(2000). 

2.1. New developments 

The experience gained from the original fortran-90 code 
allowed us to design a new version of the model, entirely 
rewritten in C++, that is much more flexible. It allows 
essential optimizations in comparison to the older model 
and a full 3-dimensional spatial grid. It is now possible to 
explicitly solve the full nuclear reaction network on a spa- 
tially resolved grid. The code can thus serve as a complete 
substitute for the conventional "leaky-box" or "weighted- 
slab" propagation models usually employed, giving many 
advantages such as the correct treatment of radioactive 
nuclei, realistic gas and source distributions etc. It also 
allows stochastic SNR sources to be included. It still con- 

4 Code WNEWTR.FOR, version 1993, as posted at 
http://spdsch.phys.lsu.edu/SPDSCH_Pages/Software_Pages/Cross_Scction_Calcs/WcbberKishSchricr.html (some minor changes have been 
made to make it compatible with GALPROP). 

5 Code YIELDX_011000.FOR, version 2000, as posted at 
http://spdsch.phys.lsu.edu/SPDSCH_Pagcs/Software_Pagcs/Cross_Scction_Calcs/SilberburgTsao.html (some minor changes have been made 
to make it compatible with GALPROP). 



tains an option to switch to the fast running cylindrically 
symmetrical model which is sufficient for many applica- 
tions such as the present one. 

In the new version, we have updated the cross-section 
code to include the latest measurements and energy de- 
pendent fitting functions. The nuclear reaction network is 
built using the Nuclear Data Sheets. The isotopic cross 
section database consists of more than 2000 points col- 
lected from sources published in 1969-1999. This includes 
a critical re-evaluation of some data and cross checks. 
The isotopic cross sections are calculated using the au- 
thor's fits to major beryllium and boron production cross 
sections p + C,N, O — > Bc,B. Other cross sections are 
calculated using Webber, Kish, & Schricr (1990) 4 and/or 
Silberberg, Tsao, & Barghouty (1998) 5 phenomenological 
approximations renormalized to the data where it exists. 
The cross sections on the He target are calculated using a 
parametrization by Ferrando et al. (1988). 

The reaction network is solved starting at the heaviest 
nuclei (i.e., Ni). The propagation equation is solved, 
computing all the resulting secondary source functions, 
and then proceeds to the nuclei with A — 1. The pro- 
cedure is repeated down to A = 1. In this way all sec- 
ondary, tertiary etc. reactions are automatically accounted 
for. This includes secondary protons, inelastically scat- 
tered primaries; their energy distribution after scattering 
is assumed to be the same as for antiprotons (eq. [2]). To 
be completely accurate for all isotopes, e.g., for some rare 
cases of /3 ± -decay, the whole loop is repeated twice. Our 
preliminary results for all cosmic ray species Z < 28 are 
given in Strong & Moskalenko (2001). 

For the calculation reported here, we use a cylindrically 
symmetrical Galactic geometry. 

2.2. Antiproton cross sections 

We calculate p production and propagation using the 
basic formalism described in Moskalenko et al. (1998). An- 
tiproton production in pp-collisions has been calculated us- 
ing the parametrization of the invariant p-production cross 
section given by Tan & Ng (1983b). This parametrization 
fits available data quite well. It gives an antiproton multi- 
plicity slightly higher than the DTUNUC code just above 
the threshold and agrees with DTUNUC results at higher 
energies (Simon et al. 1998). 

For the cross sections a™ el and a™g l we adapted 
parametrizations by Tan & Ng (1983a), Groom et al. 
(2000), and Letaw, Silberberg, & Tsao (1983). At low 
energies the total pp inelastic cross section has been calcu- 
lated using a fit from Tan & Ng (1983a). At high energies 
it is calculated as the difference between total pp and pp 
cross sections which are parametrized using Regge theory 
(Groom et al. 2000). The p absorption cross section on 
an arbitrary nuclear target has been scaled by A 2 ! 3 using 
the measured p-C,Al,Cu cross sections (Moiseev & Ormes 
1997). 
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To this we have added p annihilation and treated inelas- 
tically scattered p's as a separate "tertiary" component. 
The energy distribution after scattering is assumed to be 
(Tan & Ng 1983a) 

dN(E p ,E' p ) = 1_ 
dEp T'p 

where E'p and Ep are the total p energy before and after 
scattering correspondingly, and T p is the p kinetic energy 
before scattering. 

The p production by nuclei with Z > 2 is calculated 
using effective nuclear factors by Simon et al. (1998) and 
scaling factors similar to Gaisser & Schaefer (1992). 

In the first method we use the effective factor obtained 
from simulations of the p production with the Monte Carlo 
model DTUNUC, which appear to be more accurate than 
simple scaling. The use of this factor is consistent since the 
proton spectrum adapted in Simon et al. (1998) is close to 
our propagated spectrum above the p production thresh- 
old. For convenience, we made a fit to the ratio of the 
total p yield to the p yield from the pp-reaction (column 
3/column 2 ratio as given in Table 2 in Simon et al. 1998): 

w/vpp = 0.12 (Tp/GeV)" 1 - 67 + 1.78, (3) 

where Tp is the kinetic p energy. 

In the second method the cross section for p produc- 
tion in proton-nucleus and nucleus-nucleus interactions 
has been obtained by scaling the pp invariant cross sec- 
tion with a factor (Gaisser & Schaefer 1992) 

F u ^p X = 1.2(Aa p T l + A t a p T l )l2a p n p e \ (4) 

where Ai tt are the atomic numbers of the incident and tar- 
get nuclei, a p t, <J p i are the pAij cross sections, and a factor 
1.2 is put for consistency with Simon et al. (1998) calcu- 
lations. Production of low-energy antiprotons in pa- and 
aa-reactions is increased over that for simple scaling. In 
order to account for this effect we put Ep = _Ep + 0.06 GeV 
instead of E p when calculating the source spectra. 

Detailed results published by Simon et al. (1998) (their 
Figs. 2, 4, 7) allow us to test the approximation given by 
equation (4) and show its equivalence to DTUNUC calcu- 
lations and to equation (3). For these tests we used inter- 
stellar proton and helium spectra as published by Menn et 
al. (2000); these spectra were used by Simon et al. (1998). 
The main discrepancy is ~ 15 — 20% under-production of 
antiprotons in pp-collisions between Tp = 5 — 30 GeV com- 
pared to Simon et al. calculations (their Fig. 7), but this 
does not influence production and propagation of antipro- 
tons at lower energies. At energies below ~ 5 GeV the 
factors given by equation (3) and equation (4) are virtu- 
ally equivalent. We further use equation (3) in our calcu- 
lations. 

The pp elastic scattering is not included. At sufficiently 
high energies it is dominated by the forward peak with 
small energy transfer while at low energies the inelas- 
tic cross section (mostly annihilation) accounts for about 
70% of the total cross section (Eisenhandler et al. 1976; 
Bruckner et al. 1986). 

As we will see, the accuracy of antiproton cross sections 
is now a limiting factor; we therefore summarize the cur- 
rent status in Appendix B. 



3. SOLAR MODULATION 

Thanks to the Ulysses mission significant progress has 
been made in our understanding of the major mechanisms 
driving the modulation of CR in the hcliosphere. The prin- 
cipal factors in modulation modeling are the heliospheric 
magnetic field (HMF), the solar wind speed, the tilt of the 
heliospheric current sheet, and the diffusion tensor. Here, 
we give a short description of the model used in this pa- 
per for modulation of galactic protons and antiprotons. 
More detail on the characteristics of heliospheric modula- 
tion and modulation models can be found elsewhere (e.g., 
Potgieter 1997; Burger et al. 2000, and references therein). 

Modulation models are based on the numerical solution 
of the cosmic ray transport equation (Parker 1965): 

^M = -(v + ( Vl ,».v/ 

+ V.(K s .V/) + I(V.V)^-, (5) 

where f(r,p,t) is the CR distribution function, r is the 
position, p is the rigidity, and t is time. Terms on the 
right-hand side represent convection, gradient and curva- 
ture drifts, diffusion, and adiabatic energy losses respec- 
tively, with V the solar wind velocity. The symmetric part 
of the tensor K5 consists of diffusion coefficients parallel 
(ify) and perpendicular (K±g and K± r ) to the average 
HMF. The anti-symmetric element Ka describes gradi- 
ent and curvature drifts in the large scale HMF where the 
pitch angle averaged guiding center drift velocity for a near 
isotropic CR distribution is given by (v.d) = V x Ka&b, 
with = B/£>, and B the magnitude of the background 
magnetic field. The effective radial diffusion coefficient is 
given by K rr = K\\ cos 2 ip + K± r sin 2 ip, with ip the angle 
between the radial direction and the averaged magnetic 
field direction. We enhanced perpendicular diffusion in 
the polar direction by assuming K±g > K± r in the helio- 
spheric polar regions (see also Potgieter 2000). 

Drift models predict a clear charge-sign dependence for 
the heliospheric modulation of positively charged particles, 
e.g., CR protons and positrons, and negatively charged 
particles, e.g., electrons and antiprotons. This is due to the 
different large-scale gradient, curvature and current sheet 
drifts that charged particles experience in the HMF. For 
example, antiprotons will drift inwards primarily through 
the polar regions of the heliosphere during A < polar- 
ity cycles, when the HMF is directed towards the Sun in 
the northern hemisphere. Protons, on the other hand, 
will then drift inwards primarily through the equatorial 
regions of the heliosphere, encountering the wavy helio- 
spheric current sheet in the process. During the A > 
polarity cycles the drift directions for the two species re- 
verse, so that a clear 22-year cycle is caused (e.g., Burger 
& Potgieter 1999). 

We use a steady-state two-dimensional model that simu- 
lates the effect of a wavy current sheet (Hattingh & Burger 
1995) by using an averaged drift field with only r- and 0- 
components for the three-dimensional drift pattern in the 
region swept out by the wavy current sheet. The solar 
wind speed is 400 km/s in the equatorial plane and in- 
creases to 800 km/s in the polar regions (for details see 
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Burger, Potgieter, & Heber 2000). A modified HMF is 
used (Jokipii & Kota 1989). The outer boundary is as- 
sumed to be at 120 AU. 

The diffusion tensor is described in detail by Burger et 
al. (2000), who used the same models for turbulence (Zank, 
Matthaeus, & Smith 1996), but have adapted the coeffi- 
cients to reflect some of the results from the numerical 
simulations (Giacalonc & Jokipii 1999; Giacalone, Jokipii, 
& Kota 1999). 

The model used here gives latitudinal gradients in ex- 
cellent agreement with Ulysses observations, for both its 
value and its rigidity dependence (Burger et al. 2000). The 
model also gives realistic radial dependence for CR protons 
during consecutive solar minima, with the radial gradients 
distinctly smaller in the A > than in the A < cycles. 
Detailed fits to Pioneer and Voyager observations require 
that different diffusion coefficients must be used for con- 
secutive solar minimum periods (Potgieter 2000). 

The maximum extent of the HCS (tilt angle) has been 
used in modulation models as a proxy for solar activity 
since the 1970's, when it was realized that drifts were im- 
portant. Figure 2 shows the tilt angle 6 , which is used 
as input in our modulation calculations, for two of Hoek- 
sema's models (classic "L-model" and a newer "R-model") 
vs. time. The upper and lower curves represent his old and 
new model correspondingly (e.g., Hoeksema 1992; Zhao & 
Hocksema 1995). The difference between the models illus- 
trates the error in the tilt angles one has to consider when 
used in any modulation modeling. 



80 85 90 95 00 




80 85 90 95 00 



Years 

Fig. 2. — Heliospheric current sheet tilt angles for two Hoeksema 
models (e.g., Hoeksema 1992; Zhao & Hoeksema 1995), the classic 
"L-model" (dashes) and a newer "R-model" (solid). 

4. ANTIPROTON MEASUREMENTS 

Golden et al. (1979) and independently Bogomolovet al. 
(1979) reported the first measurements of cosmic ray an- 

6 See URL http://quake.stanford.edu/~wso/ 



tiprotons. Both experiments employed a superconducting 
magnet for measuring momentum and a gas Cherenkov 
detector to separate the background of other negatively 
charged particles from antiprotons. These were the first 
attempts to measure the flux near the 2 GeV peak ex- 
pected (Gaisser & Levy 1974) in the differential spectrum 
from secondary cosmic rays. The flux found was a factor 
of ~ 3 above that expected and caused a stir at the time. 
However, the statistics and systematic errors in both the 
measurement and the theory of that time were such that 
the excess was significant at about the 3<r level. With the 
clarity of hindsight (see below), the Golden et al. (1979) 
and Bogomolov et al. (1979) measurements were probably 
plagued by a background of negatively charged particles 
that went undetected by the Cherenkov detector. 

An annihilation technique (Buffington, Schindler, & 
Pennypacker 1981) was used to make the first alleged mea- 
surement of lower energy antiprotons (< 1 GeV). This 
technique was designed to stop the antiprotons in a visu- 
alization device. Antiprotons were identified by limiting 
their intrinsic kinetic energy (< 500 MeV) and finding 
their rest energy (~ 2 GeV) from the number of anni- 
hilation prongs and large energy release. We now know 
that this experiment must also have been plagued by a 
high background of interacting particles because the flux 
of ~ 200 MeV antiprotons found exceeds that of more re- 
cent measurements by more than an order of magnitude. 

However, these papers stimulated the realization that 
antiprotons could be produced by exotic processes such as 
the annihilation of primordial black holes or dark matter 
in the galactic halo (e.g., Carr 1985; Maki et al. 1996; Sz- 
abelski, Wdowczyk, & Wolfendale 1980; Kiraly et al. 1981; 
Silk & Srednicki 1984; Rudaz & Stecker 1988; Ellis et al. 
1988; Stecker & Tylka 1989; Jungman & Kamionkowski 
1994; Chardonnet et al. 1996). 

During the past 5 years, the BESS (Balloon-borne Ex- 
periments with a Super- conducting Solenoid) Collabora- 
tion (Yamamoto et al. 1993) has made a set of accurate 
new measurements that motivated this attempt to make a 
more accurate model calculation. This impressive payload 
now has data that extend from about 200 MeV (at the 
top of the atmosphere) to 4 GeV. Importantly, this covers 
the energy range of the expected peak in the antiproton 
distribution at 2 GeV where the flux is highest. The first 
reported mass resolved measurements (Ahlcn et al. 1988; 
Salamon et al. 1990; Yoshimura et al. 1995; Mitchell et al. 
1996) of antiprotons were made at lower energies where 
determining both the momentum per unit charge (rigid- 
ity) and velocity was easier. As can be seen from equation 
(6) both are needed to determine the mass: 

m 2 oc p 2 (/r 2 - i) , (6) 

where f3 — v/c, the velocity, and "q, the curvature, are 
the measured parameters. Rigidity, p = 1/rj, is inversely 
proportional to the curvature. Uncertainties are given by 

(Am) 2 = (Apf (A/3) 2 

m 2 p 2 (3 2 (l-j3 2 ) 2 ' 

The statistical precision in these first measurements was 
poor (Yoshimura et al. 2001). 
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In BESS the charge sign and rigidity are determined 
by track measuring drift chambers located in the uniform 
magnetic field inside the bore of their superconducting 
magnet. Extending to higher energy has required the use 
of not only the most precise time of flight measurement but 
also the addition of aerogel Cherenkov detectors for mea- 
suring velocity. The BESS team, through a series of an- 
nual balloon campaigns in Lynn Lake, Manitoba, Canada, 
has continuously improved the time resolution and in 1997 
they added an aerogel counter, extending the energy cov- 
erage beyond 2 GcV. In addition, the data handling ca- 
pacity has continuously improved to decrease dead time 
and increase precision. 

In this paper we use mainly BESS data collected in 1995 
and 1997 (Orito et al. 2000). However in one plot (§ 6, Fig. 
11), we have combined the data from those BESS flights 
with their 1998 data (Maeno et al. 2001) after correcting 
it for increased solar modulation level. The flights in 1995 
and 1997 took place near solar minimum conditions while 
the 1998 flight was just after solar modulation minimum 
(Fig. 2). We compare our predictions for different helio- 
sphcric modulation levels with data from BESS flights in 
1995-97, 1998, 1999, and 2000 (Orito et al. 2000; Maeno et 
al. 2001; Asaoka et al. 2001a) in Figures 8-10. The BESS 
flight in year 2000 provides the first accurate measurement 
of the p/p ratio for the negative polarity solar cycle. 

The statistical precision with which the flux at 2 GeV 
is known is now better than 10% (lc). The data are now 
good enough to detect solar modulation effects at this en- 
ergy at about the same magnitude. The BESS group has 
carried out an extensive calibration of their instrument 
to check for ways to reduce sources of systematic error. 
They report a systematic uncertainty in the antiproton 
flux measurements to be about 5% (Asaoka et al. 2001b). 

Data at energies above 4 GeV on the steeply falling part 
of the antiproton spectrum are more difficult to obtain. We 
have used the results reported by the MASS group (Hof 
et al. 1996; Basini et al. 1999) from a flight in 1991 based 
on the original Golden et al. (1979) payload but with an 
improved gas Cherenkov detector; we have used their most 
recent analysis (Stochaj et al. 2001). These data have er- 
ror bars that extend up and down by a factor of ~ 2. We 
can say that they are consistent with antiprotons being 
of secondary origin but are not precise enough to place 
constraints on the model presented here. 

5. NEW CALCULATIONS 

5.1. Local proton and helium spectra 

Secondary antiprotons are produced in collisions of en- 
ergetic protons and helium nuclei with interstellar gas, so 
the interstellar spectra of these nuclei is fundamental to 
calculating antiproton spectra expected at Earth. A ma- 
jor problem in determination of the LIS CR spectrum is 
the effect of heliospheric modulation that has been dis- 
cussed in detail in § 3. Inverting equation (5) is not well- 
defined making the accurate derivation of the LIS spec- 
trum a complicated task. However at energies above the 
antiproton production threshold, at 10-30 GeV/nucleon, 
the heliospheric modulation is weak. We thus try to get 
an approximate LIS spectrum using the force-field approx- 



imation (Gleeson & Axford 1968). The appropriate modu- 
lation potential (<&) has been chosen using CLIMAX neu- 
tron monitor data (Badhwar & O'Neill 1996). 

Spectra given in the literature are quoted in different 
ways, sometimes as power laws in kinetic energy, some- 
times rigidity or total energy. Also different measurements 
have different systematic uncertainties in flux. In order 
to best determine the asymptotic slope of the locally ob- 
served spectrum we have summarized all recent measure- 
ments and fitted them assuming a LIS spectrum which is 
a power-law in kinetic energy. (We note that fitting to a 
power-law in rigidity yields the spectrum which is too steep 
to agree with high energy data by Sokol and JACEE.) 

To get an idea of which part of the nucleon spectrum 
contributes most to antiproton production, we have made 
runs in which we cut the nucleon spectrum at different 
energies. Our analysis shows that ~ 97% of all antipro- 
tons below 6 GeV are produced by nucleons below ~ 200 
GeV/nucleon. The nucleons at < 20 GeV/nucleon yield 
about 1/3 of all antiprotons < 1 GeV, nucleons of < 50 
GeV/nucleon yield about 80% of all antiprotons < 2 GeV, 
and 90% of antiprotons < 4 GeV are produced by nu- 
cleons of < 100 GeV/nucleon. Therefore the spectra at 
moderate energies, up to ~ 200 GeV/nucleon, are of most 
importance. 

Tables 2 and 3 show our fits to recent measurements of 
CR protons and helium. (The power-law in kinetic energy 
has been modulated with appropriate modulation poten- 
tial and then fitted to the data.) From the parameters 
fitted to each individual set of measurements we calculate 
the weighted averages. The reduced \n (x 2 P er the degree 
of freedom) shows the quality of the fits. 

The fitted parameters are not unique since the modu- 
lation potential is not well determined, and on account 
of the systematic and statistical errors. We therefore 
use also representative high energy data, above ~ 100 
GeV/nucleon, where the modulation has no effect, for a 
cross check. In particular, our composite proton spec- 
trum, 1.60 x lO 4 ^ 75 m- 2 s- 1 sr- 1 GeV" 1 (Fig. 4, bold 
dashes), passes through Sokol (Ivancnko et al. 1993) and 
JACEE data points (Asakimori et al. 1998). The index 
agrees well with that determined by Ryan, Ormes, & Bal- 
asubrahmanyan (1972), 2.75 ± 0.03, and HEGRA below 
the knee (Arqueros et al. 2000), 2.72±°;° 2 ± 0.07. The de- 
rived LIS proton spectrum thus agrees very well with all 
the data in the range ~ 10 — 10 5 GeV. 

Our calculation of the nucleon spectrum is based on the 
solution of the diffusion equation (1), which gives the spec- 
trum in every cell of the spatial grid. Shock acceleration 
models predict that the injection CR spectrum can be ap- 
proximated by a power-law in rigidity which may steepen 
at low energies (e.g., Ellison, Slanc, & Gaensler 2001). 
Hence we adopt a power-law in rigidity, if necessary with a 
break, for the injection spectrum. The local measurements 
are conventionally used for comparison and normalization 
of the propagated spectrum. While there are some indi- 
cations that the local interstellar (LIS) nucleon spectrum 
may be not representative (Strong et al. 2000), at high 
energies where the diffusion is fast and energy losses be- 
come small it should be close to that observed near the 
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Table 2 

Parameters of the LIS proton spectrum 3, as derived in present paper. 



Instrument 


jriiLiiig iiilci vai 

-E-fcin, GeV 


1> Ul lliclllZjClLlUll, 

(m 2 s sr GeV)" 1 


± uwLi law 

index 


1V1ULI LlldLrlUIl 

potential, MV 


C It C[Llcl- 

lity, xl 


Flio-ht date 

yymmdd 


Data 

reference 


LEAP 


20-100 


— 


2.69 ±0.04 


550 


0.80 


870821 


1 


MASS2 


20-100 


(1.93 ± 0.15) x 10 4 


2.82 ±0.03 


1200 


0.53 


910923 


2 


IMAX 


20-200 


(1.08 ± 0.15) x 10 4 


2.66 ±0.04 


750 


0.26 


920716-17 


3 


CAPRICE 


20-200 


(1.55 ± 0.19) x 10 4 


2.80 ±0.03 


600 


2.54 


940808-09 


4 


AMS 


20-200 


(1.82 ± 0.21) x 10 4 


2.79 ±0.03 


550 


0.13 


9806 


5 


BESS 


20-120 


(1.61 ± 0.13) x 10 4 


2.75 ±0.03 


550 


0.05 


980729-30 


6 


Weighted average 


(1.58 ± 0.08) x 10 4 


2.76 ±0.01 










a Assuming power-law in kinetic energy LIS spectrum. 

References. — (1) Sco ct al. (1991); (2) Belotti et al. (1999); (3) Menn et al. (2000); (4) Boczio ct al. (1999); 
(5) Alcaraz ct al. (2000a); (6) Sanuki et al. (2000). 










Table 3 










Parameters of the LIS 


HELIUM SPECTRUM 3 AS DERIVED IN PRESENT PAPER. 




Instrument 


Fitting interval 
E k in, GeV/n 


Normalization, 
(m 2 s sr GeV/n)" 1 


Power-law 
index 


Modulation 
potential, MV 


Fit qua- 
lity, xl 


Flight date, 
yymmdd 


Data 
reference 


LEAP 


10-100 




2.69 ± 0.09 


550 




870821 


1 


MASS2 


10-50 


686 ± 130 


2.79 ± 0.07 


1200 


0.44 


910923 


2 


IMAX 


10-125 


600 ± 120 


2.70 ±0.07 


750 


0.69 


920716-17 


3 


CAPRICE 


10-100 


590 ± 124 


2.73 ± 0.07 


600 


0.66 


940808-09 


4 


AMS 


10-100 


653 ± 56 


2.72 ± 0.03 


550 


0.39 


9806 


5 


BESS 


10-50 


640 ± 139 


2.67 ±0.07 


550 


0.14 


980729-30 


6 


Weighted average 


641 ± 53 


2.72 ± 0.01 











a Assuming power-law in kinetic energy per nucleon LIS spectrum. 

References. — (1) Sco ct al. (1991); (2) Belotti ct al. (1999); (3) Menn ct al. (2000); (4) Boczio ct al. (1999); (5) Alcaraz ct al. (2000b); (6) 
Sanuki et al. (2000). 



Earth. This is also indicated by studies of Galactic diffuse 
emission, antiprotons, and positrons (Hunter et al. 1997; 
Moskalenko et al. 1998; Strong et al. 2000). While our 
propagation calculations yield LIS spectra which can not 
be described by a single power-law, the estimates together 
with data will serve as useful guidelines. 

The propagated helium spectrum with an appropriate 
normalization matches after modulation the data better 
than the averaged LIS spectrum from Table 3. The lat- 
ter is somewhat too low when compared with high-energy 
data by Sokol and JACEE. 

5.2. Propagation models and parameters 

To investigate the range of interstellar spectra and prop- 
agation parameters we have run a large number of models. 
Here we consider five cases (see Table 4), which differ in 
galactic propagation. The diffusive reacceleration (DR) 
model has been chosen since it has been very successful in 
the description of the propagation of CR nuclei; it repro- 
duces the sharp peak in the secondary to primary nuclei ra- 
tios in a physically understandable way without breaks in 
the diffusion coefficient and/or the injection spectrum. We 
consider also a diffusive reacceleration model with break 
in the injection spectrum (DRB) in order to match the LIS 
proton and He spectra. The minimal reacceleration plus 
convection (MRC) model uses reduced reacceleration, and 



in addition some convection is invoked to reproduce the 
B/C ratio. We compute also a plain diffusion model (PD) 
for reference, although it does not fit B/C. In the fifth 
model we use diffusion and convection (DC) with breaks 
in the diffusion coefficient and the injection index in or- 
der to construct a model capable of describing all the CR 
(particle) data simultaneously. 

In all models, the injection spectrum was chosen to 
reproduce the local CR measurements (see § 5.1). The 
source abundances of all isotopes Z < 28 are given in 
Strong & Moskalenko (2001). The propagation param- 
eters have been fixed using the B/C ratio. Because of 
the cross section fits to the main Be and B production 
channels and renormalization to the data where it exists, 
the propagation parameters show only weak dependence 
(~ 10% change in the diffusion coefficient) on the cross 
section parametrization (Silberberg et al. 1998; Webber et 
al. 1990). We thus use Webber et al. cross section code in 
our calculations. 

The halo size has been set to z h = 4 kpc, which is 
within the range zt = 3 — 7 kpc derived using the GALR- 
POP code and the combined measurements of radioactive 
isotope abundances, 10 Be, 26 Al, 36 C1, and 54 Mn (Strong 
& Moskalenko 2001). This is also in agreement with our 
newer estimate Zh — 4 — 6 kpc (Moskalenko, Mashnik, & 
Strong 2001a). Note that the exact value of Zh is unimpor- 
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Table 4 
Propagation parameter sets. 



Diffusion cocfficicnt a Rcaccclcration'/Convcction 



iviociei 


tnjection index, r y c 


Do , cm 2 s ^ 


Index, 8 


va/w^/'^, km s ^ 


dV/dz, km s ^ kpc ^ 


Diffusive 

Rcaccclcration (DR.) 


O A Q 

z.4o 


0. 1U X 1U 


U.oo 






Diffusive 
Rcacceleration 
with Break (DRB) 


1.93/2.43 d 


6.10 x 10 28 


0.33 


30 




Minimal Reacceleration 
plus Convection (MRC) 


2.43 


4.30 x 10 28 


0.33 


17 


10 


Plain Diffusion (PD) 


2.16 


3.10 x 10 28 


0.60 






Diffusion plus 
Convection (DC) 


2.46/2.16° 


2.50 x 10 28 


0/0.60 a 




10 



Note. — Adopted halo size Zh = 4 kpc. 

a p0 = 4 GV, index 5 is shown below/above po- 

b «A is the Alfven speed, and w is denned as the ratio of MHD wave energy density to magnetic field energy density. 



c For a power-law in rigidity, <x p 1 . 
d Index below/above rigidity 10 GV. 
e Index below/above rigidity 20 GV. 



tant for antiproton calculations provided that the propa- 
gation parameters are tuned to match the B/C ratio. 

Our results are plotted in Figures 3-7. The upper curve 
is always the LIS spectrum, except in the B/C ratio plot 
(Fig. 3) where the lower curves are the LIS ratio. The 
modulation on these plots has been done using the force- 
field approximation. For the B/C ratio (Fig. 3) in the DC 
model the result of modulation in the drift model is also 
shown though the difference (force-field approximation vs. 
drift model) is small. 

The "tertiary" antiprotons (inelastically scattered sec- 
ondaries), significant at the lowest energies, are important 
in interstellar space, but make no difference when com- 
pared with measurements in the heliosphere (Fig. 7). 

The model with the diffusive reacceleration (DR) re- 
produces the sharp peak in secondary to primary nuclei 
ratios (Strong & Moskalcnko 1998, 2001). However, this 
model produces a bump in proton and He spectra at <~ 2 
GeV/nucleon that is not observed. 7 This bump can be 
removed by choosing an injection spectrum that is flat- 
ter at low energies (see, e.g., Jones et al. 2001). There are 
however some problems with secondaries such as positrons 
and antiprotons which are difficult to manage. A similar 
bump appears in the positron spectrum at ~ 1 GeV, and 
the model underproduces antiprotons at 2 GeV by more 
than 30%. 

In order to fit the proton/Helium spectra in a diffu- 
sive reacceleration model we need to introduce a break in 
the injection spectrum, with a smaller index below 10 GV 
(model DRB, not shown in the plot for clarity). This break 
in the injection spectrum of primaries results in smaller 
bump in the secondary posiron spectrum (dash-dots in 
Fig. 6), though the predicted spectrum still shows a sig- 
nificant excess of a factor at least two over the measure- 
ments. However, it produces no effect on secondary an- 



tiprotons because of their higher production threshold en- 
ergy. Taken together these cases provide evidence against 
"strong" reacceleration 8 in the ISM. 




I i i— lJ i i— iJ i i— iJ i i— iJ i i— J 

0.01 0.1 1 10 100 1000 

Kinetic energy, GeV/nucleon 

Fig. 3. — B/C ratio as calculated in DC, DR, and PD models 
for Zh = 4 kpc. Lower curves — interstellar, upper — modulated 
($ = 450 MV). The line coding is: solid lines - DC model, dashes 
— DR, dots - PD. The thick solid line (DC model) corresponds to 
the drift model calculations of hcliosphcric modulation for A > 
and tilt angle 5°. Data below 200 MeV/nucleon: ACE (Davis ct 
al. 2000), Ulysses (DuVernois, Simpson, & Thayer 1996), Voyager 
(Lukasiak, McDonald, & Webber 1999); high energy data: HEAO-3 
(Engclmann ct al. 1990), for other references see Stephens & Stre- 
itmatter (1998). 



7 A similar bump is produced also in the electron spectrum at 1 GeV. 

8 We define that the reacceleration is "strong" if the model is able to match the B/C ratio without invoking other mechanisms such as 
convection and/or breaks in the diffusion coefficient. 
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Kinetic energy, GeV 
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Fig. 4. — Calculated proton interstellar spectrum (LIS) and modulated spectrum (force field, <3? = 550 MV). The lines are coded as in 
Figure 3. Bold dashed line shows our composite spectrum as derived from measurements 1.60 X 10 4 E^ n 7 '' m -2 s _1 sr _1 GeV -1 (not shown 
between 5-500 GeV for clarity). Data: IMAX (Menn et al. 2000), CAPRICE (Boezio ct al. 1999), AMs"(Alcaraz et al. 2000a), BESS (Sanuki 
et al. 2000), Sokol (Ivanenko et al. 1993), and JACEE (Asakimori et al. 1998). 




10 100 
Kinetic energy, GeV/nucleon 



1000 



Fig. 5. — Calculated He interstellar spectrum (LIS) and modulated spectrum (force field, <3? = 550 MV). The lines are coded as in Figure 
3. Data: IMAX (Menn ct al. 2000), CAPRICE (Boezio ct al. 1999), AMS (Alcaraz et al. 2000b), BESS (Sanuki ct al. 2000), Sokol (Ivanenko 
et al. 1993), and JACEE (Asakimori et al. 1998). 
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Fig. 6. — Calculated positron flux: upper curves LIS, lower mod- 
ulated (force field, <t> = 600 MV). The lines are coded as in Figure 
3. Positron flux as calculated in DRB model is shown by dash-dots. 
Data: HEAT94 (Barwick ct al. 1998), CAPRICE94 (Boczio et al. 
2000). 
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* CAPRICE98 
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ANTIPROTONS 

Phi = 550 MV 




0.001 



Kinetic energy, GeV 

Fig. 7. — Calculated antiproton spectrum, upper curves - LIS, 
modulation was made with $ = 550 MV (force field, lower curves). 
The lines are coded as in Figure 3. Data: BESS (Orito ct al. 2000), 
MASS91 (Stochaj ct al. 2001), CAPPJCE98 (Boezio et al. 2001). 
The top curves are the total secondary p's. The two lowest curves 
marked "tertiary" show separately the LIS spectrum and modulated 
"tertiary" component in DC model. 



Another model combining reduced reacceleration and 
convection (MRC) also produces too few antiprotons. 

Using a plain diffusion model (PD) we can get good 
agreement with B/C above few GeV/nucleon, with nu- 
cleoli spectra and positrons, but this model overproduces 
antiprotons at 2 GeV by ~ 20% and contradicts the sec- 
ondary/primary nuclei ratio (B/C) below 1 GeV/nucleon. 

A diffusion model with convection (DC) is our "best 
fitting model". It reproduces all the particle data "on 
average" , although it has still some problem with repro- 
duction of the sharp peak in the B/C ratio. In this model 
a flattening of the diffusion coefficient (S = 0) below 4 
GV is required to match the B/C ratio at low energies. A 
possible physical origin for the behaviour of the diffusion 
coefficient is discussed in § 5.3. 

To better match primaries (jp, He) in the DC model, we 
introduced a steeper injection spectrum below 20 GV to 
compensate the break in the diffusion coefficient, which 
tends to flatten the spectrum at low energies. Such a 
steepening in the injection spectrum, however, has almost 
no effect on secondaries (p, e + ). The existence of a sharp 
upturn below few GeV/nucleon is in fact predicted from 
SNR shock acceleration theory (Ellison, Slane, & Gaensler 
2001); this is a transition region between thermal and non- 
thermal particle populations in the shock. Our model does 
not require a large break (0.3 in index is enough, see Table 
4). More discussion is given in § 5.3 and § 7. 

The drift model has been used for heliospheric modula- 
tion of all spectra obtained in DC model (p, He, B/C, e+, 



p). The same modulation parameter set has been applied 
to all the species and the agreement is good. The B/C ra- 
tio (Fig. 3) is virtually insensitive to the modulation level, 
and it remains the same for a wide variety of tilt angles. 

5.3. Discussion 

Our "best-fit" model (DC) with a constant diffusion co- 
efficient below 4 GV suggests some change in the propaga- 
tion mode. Propagation and scattering of the high energy 
particles on magnetic turbulence is described by diffusion 
with the diffusion coefficient increasing with energy. The 
growth of the diffusion coefficient depends on the adopted 
spectrum of magnetic turbulence and is typically in the 
range 0.3 — 0.6. However, at low energies particles could 
propagate following the magnetic field lines rather than 
scatter on magnetic turbulence. This change in the prop- 
agation mode may affect the diffusion coefficient making 
it less dependent on energy. Since the magnetic field lines 
are essentially tangled such a process can still behave like 
diffusion (Berezinskii et al. 1990). 

One more unknown variable is the CR spectrum in the 
distant regions of the Galaxy. The LIS nucleon spectrum 
is studied quite well by direct measurements at high en- 
ergies where solar modulation effects are minimal. Mean- 
while the ambient CR proton spectrum on the large scale 
remains unknown. The most direct test is provided by dif- 
fuse 7-rays. However there is the well known puzzle of the 
GeV excess in the EGRET diffuse 7-ray spectrum (Hunter 
et al. 1997) which makes a direct interpretation in terms of 
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protons problematic, and possibly inverse-Compton emis- 
sion from a hard electron spectrum is responsible (Strong 
et al. 2000); for this reason we do not consider 7-rays fur- 
ther here. 

An interesting independent possibility to test the inter- 
stellar nucleon spectrum at different energies is provided 
by antiproton and positron data taken together. While the 
spectra of antiprotons and positrons at above ~ 1 GeV al- 
low conclusions on the nucleon spectrum above ~ 10 GeV 
(Moskalenko et al. 1998), the positron spectrum below 1 
GeV is sensitive to the nucleon spectrum below the an- 
tiproton production threshold ~ 10 GeV. The problem 
here is that the propagation parameters are not known 
accurately enough and this affects the predicted spectrum 
of positrons at low energies. 

6. VARIATIONS OF PROTON AND ANTIPROTON SPECTRA 
OVER THE SOLAR CYCLES 

We use the DC model to calculate the LIS spectra of 
protons and antiprotons and then use the drift model to 
determine their modulated spectra and ratio over the solar 
cycles with positive (A > 0) and negative (A < 0) polar- 
ity (Figs. 8-10). The variations shown depend on the tilt 
angle. When combined with Figure 2 this allows us to es- 
timate the near-Earth spectra for arbitrary epochs in the 
past as well as make some predictions for the future. It 
may be also used to test the theory of heliospheric modu- 
lation. 

The p/p ratio is supposed to be more accurately mea- 
sured than the flux. Interestingly, the p/p ratio ap- 
pears largely insensitive to the modulation level during 
the A > cycle (Fig. 10, left). This agrees within the er- 
ror bars with BESS p/p measurements made in 1993, 1995, 
1997, 1998 (Maeno et al. 2001), and 1999 (Asaoka et al. 
2001a), and thus allows us to constrain the proton spec- 
trum at low energies, in the heliosphere above ~ 0.1 GeV 
and the LIS spectrum above ~ 0.7 GeV. In the new cycle 
(A < 0), on the contrary, this ratio is predicted to vary 
by over an order of magnitude (Fig. 10, right), which will 
allow to us validate the calculations of heliospheric modu- 
lation including charge sign effects. Our predictions for tilt 
angles 55° — 65° (indicated by the arrow) agree well with 
the data from the BESS-2000 flight (Asaoka et al. 2001a). 
The HEAT-2000 data (Musser et al. 2001) were obtained 
during the Spring 2000 flight at around the solar maxi- 
mum, but the solar modulation is weak at high energies. 
Shown are also two data points obtained at the moderate 
level of solar activity in 1984-1985 and 1986-1989 in the 
previous A < cycle (Bogomolov et al. 1987, 1990, and 
Bogomolov, private communication). 

The proton flux data obtained during BESS-1999, 2000 
flights (Asaoka et al. 2001a) agree qualitatively with varia- 
tions over the solar cycle. However, BESS-1999 data show 
a flux which is somewhat larger while BESS-2000 data 
show the flux lower than predicted. Both measurements 
have been made during maximum solar activity and near 
the solar magnetic field reversal and may reflect the dy- 
namic effects taking place in the heliosphere. 

Now we use the predicted variations over the solar cycle 
to combine BESS data collected in 1995, 1997 during the 



solar minimum with BESS data collected in 1998, when 
the solar activity was moderate (the tilt angle about 25°). 
The correction was calculated using the energy-dependent 
factor, F(E) 5 /25 = ^5/^25, the ratio of modulated an- 
tiproton spectra with tilt angles 5° and 25°. The data 
of 1998 and their error bars have been multiplied by this 
factor. The spectrum obtained in such a way has been 
further combined with BESS data collected in 1995, 1997. 

The combined spectrum agrees very well with calcula- 
tions for tilt angle 5° and A > (Fig. 11). It is, however, 
not sensitive to the exact value of the tilt angle at the 
epoch of BESS-98 flight. It looks very similar for tilt an- 
gles 15° and 35°. The reason is that the error bars of 
the data collected in 1998 are large, so that their relative 
weight is small. 

7. CONCLUSIONS 

In this paper, we have studied five basic propagation 
models. Table 5 summarizes the results of this study. 
First, we applied a propagation model with reacceleration. 
It reproduces nuclear secondary/primary ratios rather well 
but has some difficulties with reproduction of primary pro- 
ton and He spectra, and, more important with secondary 
positrons and antiprotons. When the primary injection 
spectra arc adjusted to fit the observed modulated p and 
He spectra, the secondary positron excess is reduced but 
still significant, while the antiproton deficit is unaffected. 
Another model with reduced reacceleration strength and 
convection also produces too few antiprotons. A plain dif- 
fusion model (no reacceleration, no convection) reproduces 
the high energy part of nuclear secondary /primary ratios, 
protons, and positrons, but it has a problem with the low 
energy part of nuclear secondary/primary ratios and over- 
estimates the antiproton flux. A model with convection 
and flattening of the diffusion coefficient at low energies 
and a break in the injection spectrum is our "best-fit" 
model. It reproduces well all particle data "on average". 
The spectrum of diffuse 7-rays calculated for this model, 
cannot explain the EGRET GeV excess (Hunter et al. 
1997), but this could originate in other ways. 

During the last decade there have been a number of 
space and balloon experiments with improved sensivity 
and statistics. They impose stricter constraints on the 
CR propagation models. It has become clear that cur- 
rently there is no simple model that is able to simultane- 
ously reproduce all data related to CR origin and prop- 
agation. This conclusion is mainly the result of the in- 
creased precision of the CR experimental data, but also 
the improved reliability of the calculations of interstellar 
and heliospheric propagation. 

What could be the origin of this failure (to find a simple 
model), apart from the propagation models? Concern- 
ing the accuracy of the experimental data, the spectra of 
protons and helium arc measured almost simultaneously 
and quite precisely by BESS and AMS and the agreement 
is impressive. They also agree with earlier experiments, 
such as LEAP, IMAX, CAPRICE, within the error bars. 
The most accurate measurements of nuclei at low energies 
are made by Voyager, Ulysses, and ACE, and the agree- 
ment is good. At higher energies the data obtained by 
HEAO-3 are the most accurate and generally agree with 
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Fig. 8. — Calculated proton LIS and modulated spectra for the two magnetic polarity dependent modulation epochs, A > (top) and A < 
(bottom). Tilt angle from top to bottom: 5°, 15°, 25°, 35°, 45°, 55°, 65°, 75°. The tilt angle corresponding to BESS-98 data is ~ 5° - 15° 
(A > 0) depending on the coronal field model. For discussion of BESS-99 and BESS-00 data see text. Data: BESS-98 (Sanuki et al. 2000), 
BESS-99, BESS-00 (Asaoka et al. 2001a). 
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Fig. 9. — Calculated antiproton LIS and modulated spectra for the two magnetic polarity dependent modulation epochs, A > (top) and 
A < (bottom). Tilt angle from top to bottom: 5° , 15° , 25° , 35° , 45° , 55° , 65°, 75°. The tilt angle corresponding to BESS data is ~ 5° - 15° 
(A > 0) depending on the coronal field model. Data references as in Figure 7, BESS-99, BESS-00 (Asaoka et al. 2001a). 
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Fig. 10. — Calculated p/p ratio in the interstellar medium (LIS) and modulated for the two magnetic polarity dependent modulation epochs, 
A > and A < 0. Top: Tilt angle from top to bottom (solid lines): 5°, 15°, 25°, 35°; from top to bottom (dotted lines): 75° , 65° , 55° , 45° . 
(Lines almost coincide for the tilt angles 35°, 45° and 25°, 55°, and very close for 15°, 65°.) Data: BESS (Orito et al. 2000; Maeno et al. 
2001; Asaoka et al. 2001a), MASS91 (Hof et al. 1996; Stochaj et al. 2001), and CAPRICE98 (Bergstrom et al. 2000). Bottom: Tilt angle 
from bottom to top (solid lines) : 5° , 15° , 25° , 35° , 45° ; from bottom to top (dotted lines) : 75° , 65° , 55° . (Lines are very close for the tilt angles 
35°, 65° and 45°, 55°.) The BESS-00 data agree well with our predictions for the tilt angles 55° - 65° (indicated by the arrow). Data: IPTI 
(Bogomolov et al. 1987, 1990), HEAT-00 (Musser et al. 2001), BESS-00 (Asaoka et al. 2001a). 
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Table 5 

Summary of model predictions. 



Primaries Secondaries 



Model 


B/C 


Protons 


He 


Antiprotons 


Positrons 


DR a 


Good 


LE bump/Good 


Fair / Good 


Too few 


LE bump 


DRB 


Good 


Good 


Good 


Too few 


LE excess 


DR a 


Good 


LE bump/Good 


Fair / Good 


Too few 


LE bump 


MRC a 


Fair 


LE bump/Good 


Fair / Good 


Too few 


LE bump 


PD 


Too large at LE 


Good 


Good 


Too many 


Good 


DC 


Fair 


Good 


Good 


Good 


Good 



a The reproduction of spectra of primaries can be improved by choice of the injection spectrum. 



earlier measurements. Electron flux measurements made 
by HEAT, CAPRICE, and at Sanriku balloon facility all 
agree. Positron data, though with large error bars, are 
in agreement as well. The possibility that BESS antipro- 
ton data have systematic errors as large as +20% looks 
improbable. They have carried out calibrations which 
demonstrate that the error must be < 5% (Asaoka et al. 
2001b). However, additional measurements are desirable. 




0.1 1 10 

Kinetic energy, GeV 



Fig. 11. — BESS 1995-97 data and the data combined with "cor- 
rected" BESS 1998 measurements. Shown also are the calculated 
antiproton interstellar (LIS) and modulated spectrum for a tilt angle 
5° (A > 0). Data references as in Figure 7. 

Nuclear cross section errors are one of the main con- 
cerns. Fitting (matching) the measured B/C ratio is a 
standard procedure to derive the propagation parameters. 
The calculated ratio, in turn, depends on many cross sec- 
tions, such as the total interaction and isotope production 
cross sections. The latter appear to have quite large un- 
certainties, typically about 20%, and sometimes they can 
even be wrong by an order of magnitude (see, e.g., Mos- 
kalenko et al. 2001a). This is reflected in the value of 
the diffusion coefficient, the Alfven velocity (parameter- 



izing reacceleration) and/or convection velocity, and thus 
influences the calculated spectra of CR species. Our cross 
section calculations make use of the fits to the cross sec- 
tions p + C,N, O — ► Bc,B, that produce most of the Be 
and B. Other channels are calculated using the Webber et 
al. or Silberberg et al. cross section codes renormalized to 
data where it exists. We thus can rule out a possibility of 
large errors in the calculated B/C ratio. 

Solar modulation for antiprotons is different from that of 
protons, the effect known as charge sign dependence. Over 
the last years Ulysses made its measurements at different 
heliolatitudes so we know more about the solar magnetic 
field configuration and the solar wind velocity distribution. 
On the other hand, the two Voyagers, the most distant 
spacecraft, provide us with measurements of particle fluxes 
made close to the heliospheric boundary. There is thus a 
small chance of a serious error, e.g., that the antiproton 
flux during the last solar minimum was modulated much 
more weakly than estimated. 

The underproduction of antiprotons in the reaccelera- 
tion model may be connected with a contribution of pri- 
mary antiprotons, but this suggestion conflicts with other 
CR data. A strong primary antiproton signal would be ac- 
companied by the positron signal and an excess in 7-rays 
(EGRET GeV excess ?). However, the positron flux cal- 
culated in the DR model already shows an excess rather 
than deficit. Besides, only few SUSY dark matter can- 
didates are able to produce a signal large enough to be 
detected and even in this case primary antiprotons con- 
tribute mostly to low energies (see discussions in Bottino 
et al. 1998; Bergstrom et al. 1999). 

Considering reacceleration models, there is the possibil- 
ity that the injection spectra of protons and nuclei might 
be different. For instance, spectra of nuclei Z > 2 are well 
reproduced by a reacceleration model with a power-law (in 
rigidity) injection spectrum, while protons, helium, and 
electrons require some spectral flattening at low energies 
to avoid developing a bump in the spectrum. This is in 
agreement with conclusions made by other authors (Jones 
et al. 2001). The difficulty is to get the right positron 
and antiproton fluxes in the reacceleration model. On the 
other hand, non-reacceleration models require some spec- 
tral steepening at low energies (compensation for flatten- 
ing of the diffusion coefficient) to get good agreement with 
the spectra of primaries. Such steepening may be a con- 
sequence of SNR shock acceleration (Ellison et al. 2001). 
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If we assume that the problem is in propagation mod- 
els, we have shown that it is possible to construct a model 
(DC) that fits all these data, by postulating a significant 
flattening of the diffusion coefficient below 4 GV together 
with convection (and possibly with some reacceleration) 
and a break in the injection spectrum. The break in the 
diffusion coefficient is reminiscent of the standard proce- 
dure in "leaky-box" models where the escape time is set 
to a constant below a few GeV. This has always appeared 
a completely ad hoc device without physical justification, 
but the present analysis suggest it may be forced on us 
by the data. Therefore, possibilities for its physical origin 
should be studied. 

Where do we go from here? New measurements of CR 
species that cover the range 0.5 — 1000 GeV/nucleon are 
necessary to distinguish between reacceleration and non- 
reacceleration models. The new experiment PAMELA 
scheduled for launch in 2002 should improve accuracy of 
positron (0.1-200 GeV) and electron (0.1-300 GeV) mea- 
surements while allowing for simultaneous measurements 
of p and nuclei from H to C in the energy range 0.1-200 
GeV/nuclcon (The PAMELA Coll. 1999). The puzzling 
excess in the diffuse 7-rays above 1 GeV in EGRET data 
should be confirmed. Here we await the launch of GLAST 
with 30 times better sensivity. Its improved angular res- 
olution will be helpful in removing the unresolved source 
component contribution to the diffuse flux. A possibil- 
ity that the nucleon spectrum in the Galaxy is harder or 
softer than we measure locally should be explored by mak- 
ing more accurate measurements of antiproton fluxes at 
high energies (Moskalenko et al. 1998). This may help to 
produce more antiprotons in the DR model and simulta- 
neously explain the GeV excess in 7-rays. The interpreta- 
tion of B/C and other secondary /primary ratios depends 
on the assumption that all the secondaries are produced 
during propagation; if instead some secondary production 
occurs in the sources the range of possible models widens. 
Whether such effects could provide an alternative explana- 



tion of the observed effects should be the subject of future 
study. The effect of the local distribution of interstellar 
matter is worth investigating and we are going to do this 
in future by developing and incorporating in our code a 
3D model of the gas distribution in the Galaxy. The ef- 
fect should however only important for radioactive and 
K-capture isotopes and electrons/positrons, i.e. for those 
particles whose half-life is less or comparable to the cos- 
mic ray life-time in the Galaxy or whose large energy losses 
confine them to a short distance from the source. On the 
other hand for protons, with their 30 mb total cross section 
and negligible energy losses, only the large-scale gas distri- 
bution is important. The same is true for antiprotons, as 
their total interaction cross section at 1 GeV is only about 
two times larger and approaches the proton cross section 
at higher energies (Figure B12). A thorough and system- 
atic measurement of the nuclear cross sections over a wide 
energy range is extremely important. Accurate measure- 
ments of the cross sections would allow more information 
to be extracted from the now rather precise CR measure- 
ments. New measurements of the antiproton and charged 
and neutral pion production cross sections are desirable. 
The widely used parametrization of antiproton production 
and interaction cross sections by Tan & Ng (1983a, b) is 
twenty years old. Parametrization of pion production is 
even older. Given the new experimental techniques now 
available such measurements should be straightforward. 
At the same time the sensitivity of the whole analysis to 
the modulation models has to be investigated further. 
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APPENDIX 

THE INTERSTELLAR GAS DENSITY DISTRIBUTION IN THE CYLINDRICALLY SYMMETRICAL MODEL 

The H2 number density in mols. cm~ 3 is calculated from 

n K2 {R,z) = 3.24 x 10 - 22 Ae (i?) e~ ln2 ^-^/ z l 7 (Al) 

where cq(R) (K km s _1 kpc -1 ) is the CO volume emissivity, z (R) and Zh(R) are the height scale and width defined by 
a table (Bronfman et al. 1988), and X = raH 2 / e co = 1-9 x 10 20 mols. cm~ 2 /(K km s -1 ) is the conversion factor (Strong 
& Mattox 1996). 

The H I (atom cm -3 ) relative distribution is taken from Gordon & Burton (1976), but renormalized to agree with 
Dickey & Lockman (1990), since they give their best model for the ^-distribution in the range R — 4 — 8 kpc and state 
the total integral perpendicular to the plane is 6.2 x 10 20 cm -2 : 

1 f E J= i2^ e ~ ln2 z2/z * + 4e~ |z|/z3 , i?<8kpc 

n H i(R,z)= Y(R) I interpolated, 8< J R<10kpc (A2) 

nGB { n DL cxp(-z 2 e- - 22i 7z 2 ), R > 10 kpc 

Here Y(R) is the distribution from Gordon & Burton (1976) (R < 16 kpc), jigb = 0.33 cm -3 and tidl = 0.57 cm -3 are 
the disk densities in the range 4 < R < 8 kpc in models by Gordon & Burton (1976) and Dickey & Lockman (1990), 
correspondingly. The z-dependence is calculated using the approximation by Dickey & Lockman (1990) for R < 8 kpc, 
using the approximation by Cox et al. (1986) for R > 10 kpc, and interpolated in between, and the parameter values are 
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A 1 = 0.395, A 2 = 0.107, A 3 = 0.064, z x = 0.106, z 2 = 0.265, z 3 = 0.403, z A = 0.0523. For R > 16 kpc an exponential tail 
is assumed with scale length 3 kpc. 

The ionized component H II (atom cm -3 ) is calculated using a cylindrically symmetrical model (Cordes et al. 1991): 

n H u(R,z) = £ n t e-^ h ^ R -^ 2 / a l (A3) 

i=l,2 

where n x = 0.025, n 2 = 0.200, hi — 1 kpc, h 2 = 0.15 kpc, R x = 0, i? 2 = 4 kpc, a! = 20 kpc, a 2 = 2 kpc. 

PROTON AND ANTIPROTON CROSS SECTIONS 

The energy and momentum units are GeV and GeV/c correspondingly. All cross section given are in mbarn and plotted 
in Figure B12 (except for the production cross section). An asterisk marks the center-of-mass system (CMS) variables. 
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Fig. B12. — Proton and antiproton cross sections. 



The inclusive antiproton production cross section (mb GeV c 3 ) in pp-reaction is given by (Tan & Ng 1983b): 

E P-TT = S{x t )f(x r ) exp{-A(x r )p t + B(x r )Pt}, 
dpp 



(Bl) 



/ = 1.05 x 10~ 4 cxp(-lO.lav) 0(0.5 - x r ) + 3.15(1 - x r ) 7 9 , 

A = 0.465 cxp(-0.037x r ) + 2.31 exp(0.014x r ), 

B = 0.0302 cxp{-3.19(av + 0.399)}(x r + 0.399) 8 39 , 

x r = E;/E* p max , 

where 9(x) is the Heaviside step function (9\x > 0] = 1, otherwise =0), pt is the transverse momentum of the antiproton, 
Ep is the total CMS energy of the antiproton, Ep max is the maximal value of Ep for the given inclusive reaction, and 

S(xt) is the low energy correction, S = 1 at s 1 / 2 > 10 GeV. At s 1 / 2 < 10 GeV it is given by: 

S- 1 = 1 - cxp {- cxp [c(x t )Q - d(x t )] (l - cxp [-a{x t )Q b ^ ) } , 

a = 0.306 exp(-0.12a; t ), 
0.0552 cxp(2.72a; t ), 

0.758 - 0.68a; t + 1.54a; 2 , (B2) 
0.594 exp(2.87x t ), 



b- 
c ■ 
d- 



3 1 ' 2 - AM, 



pi 



ti* jrp* max 

t pip i 

s = 2M p (Ep + M p ) = inv is the square of the total CMS energy of colliding particles, M p is the proton rest mass, T~ is 
the kinetic CMS energy of antiproton, and Tp max is the maximal value of Tp for the given inclusive reaction. 
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Proton-proton inelastic cross section (Tan & Ng 1983a): 

T p < 0.3; 

Ariel _ on o h i n no7irr i n m n^amw ) < 1 , n on . . 1 n-3r^-C\ _1 



a™ e( = 32.2 [1 + 0.0273*7 + O.O1C/ 2 0(J7)] <| (l + 2.62 x IQ- 3 T- C ) , 0.3 < T p < 3; 
f/ = ln(£ p /200), 



T P > 3; 



C* = 17.9+ 13.8 lnT p + 4.41 In 2 T p , 



where #(£/) is the Heaviside step function, £7 p and T p are the total and kinetic energy of proton, correspondingly. 
Nucleus-proton inelastic cross section (Letaw et al. 1983): 



a p X = 45A U 7 5(T p ) [1 + 0.016 sin(5.3 - 2.63 In A)] 



1 - 0.62e- T ?/°- 2 sin 



10.9 



Tp ^ 3, 



1, T p > 3; 



5 = 



1 + 0.75exp{-T p /0.075}, for beryllium; 
1, otherwise; 



where A is the atomic number. In the case of pHe inelastic cross section equation (B3) is known to be not very accurate, 
while the p— 4 He cross section is the most important after the pp cross section. We thus made our own fit to the data in 
the range 0.02 - 50 GeV: 

4 n He = HI {1 - cxp[-3.84 (T p - 0.1)] (1 - sin[9.721og a319 (10 3 T p ) - 4.14])} , (B3) 

where T p > 0.01 GeV. 

The antiproton-proton annihilation cross section at low energies is taken from Tan & Ng (1983a). At high energies it is 
calculated as the difference between the total pp and pp cross sections which are parametrized using Regge theory (Groom 
et al. 2000): 

661 (1 + 0.0115 T p - 774 - 0.948 T? 0151 ) , T p < 10; , , 

2 x 35.43s-°- 56 , T p > 10; { } 

where s is the square of the total CMS energy of colliding particles as defined in equation (B2). Annihilation is of minor 
importance at high energies. 

Total inelastic antiproton-proton cross section (Tan & Ng 1983a): 

tot J 24.7 (l + 0.584 Tf - 115 + 0.856 Tp - 566 ) , T p < 14; 

° pp <y a ,nel +(J a^ T p > 11. ^> 

At very low energies, below ~ 0.01 GeV, we put Op p = (7*°,*. 

Non- annihilation inelastic pp cross section is calculated as Op° n = a l p p — o°p p . The antiproton- nucleus non- annihilation 
inelastic cross section is taken equal to proton- nucleus inelastic cross section, a™^ = a p n ^ . 

The parametrization of the p total inelastic cross section on an arbitrary nuclear target has been obtained in Moiseev 
& Ormcs (1997) using a parametrization by Kuzichev, Lcpikhin, & Smirnitsky (1994) and has been tested against the 
data available on p cross sections on C, Al, and Cu targets: 

afl = A 2/3 [48.2 + 19 T p a - 55 + (0.1 - 0.18 T p 12 )Z + 0.0012 Tf^Z 2 ] , (B6) 

where Z is the nucleus charge. The second term in square brackets has been modified from its original form 19 (Tp — 
0.02)~ 55 to better match the slope below ~ 100 MeV. This modification does not affect agreement with data at higher 
energies. Another modification, in the case of 4 He, is to use A = 3.3 to scale the value down slightly to be consistent with 
cross section data both at low and at high energies (Allaby et al. 1971; Balcstra et al. 1985). 

The pA annihilation cross section is calculated from Op A = a p °X — ^"ai where a™ A n = (J p n A el . This is quite accurate 
provided that the He abundance is about 10% that of H. 
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